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ABSTRACT 

A fraction of the mid-J (J= 14-13 to J— 24-23) CO emission detected by 
the Herschel /PACS observations of embedded young stellar objects (YSOs) has 
been attributed to the UV-heated outflow cavity walls. We have applied our 
newly developed self-consistent models of Photon-Dominated Region (PDR) and 
non-local thermal equilibrium line Radiative transfer In general Grid (RIG) to 
the Herschel FIR observations of 27 low mass YSOs and one intermediate mass 
YSO, NGC7129-FIRS2. When the contribution of the hot component (traced 
by transitions of J > 24) is removed, the rotational temperature of the warm 
component is nearly constant with ~ 250 K. This can be reproduced by the 
outflow cavity wall (n > 10 6 cm -3 , logG 0 /n > —4.5, log G 0 > 3, T gas > 300 K, 
and X(CO)> ICE 5 ) heated by a UV radiation held with a black body temperature 
of 15,000 K or 10,000 K. However, a shock model combined with an internal PDR 
will be required to determine the quantitative contribution of a PDR relative to 
a shock to the mid-J CO emission. 


1. Introduction 

Embedded young stellar objects (YSOs) are associated with energetic phenomena: jets, 
outflows, and high energy photons emitted by accretion shocks on the surface of protostars 
and disks. These phenomena determine the physical conditions of the surrounding material, 
in particular, in close proximity to the central object. However, it is difficult to directly 
detect emission from the warm/hot gas and dust closest to the forming star because of the 
thick enshrouding envelope. 

In this regard, far infrared (FIR) spectroscopy can be a powerful tool for the study 
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of embedded YSOs because the energetic photons produced by accretion are absorbed 
and re-emitted in this wavelength re gime. FIR spec t 


low-mass embedded protostars (e.g., 


Giannini et ah 


jtroscc 

pic observations of twentv-eight 

2001; 

Nisini et al. 

2002; 

van Dishoeck 


20Q4r) were carried out for the first time with the Long Wavelength Spectrometer (LWS; 


Clegg et ah 


1996) aboard the Infrared Space Observatory (ISO). These observations 


discovered widespread CO emission arising f rom a rotational s tates from J= 14 to J— 29 


(for NGC1333 IRAS 4, 


Giannini et al 


2001 


Maret et al. 


2002). A Large Velocity Gradient 


(LVG) analysis of these observations suggests that the FIR CO emission is radiated from 

i n5 „, i n6 


(e-g., 


s with the term 

)erature of a few hundret 

:s to ~1000 K and a 

Giannini et al. 

2001; 

Nisini et al. 

2002; 

van Dishoeck 

to 

o 

o 

^ _ 


2004). A non-dissociative shock 


is the leading candidat e for the heating mechanism pro ducing the high-J CO emission 


(e-g-, 


Nisini et al. 


5 ; 


2002b however 


Ceccarelli et al. 


(120021 1 proposed that super-heated a disk 


surface layer that is exposed to ultraviolet photons could account for the FIR CO emission 
in Elias 29. 


The Herschel Space Observatory (jPilbratt et al. 


2010) provides higher spatial resolution 


and sensitivity when compared to ISC 

. Furthermore the Photodetector Array Camera and 

Spectrometer (PACS, 

Pogli 

tsch et al. 

2010 

) covers the CO rotational lines from J= 14-13 

to J= 49-48 ( 

Herczeg et al. 

2012 

) and provided a wealth of new observations of energetic 


gas in protostars. A Herschel open time key program, “Herschel Orion Protostar Survey 


(HOPS)” observed 22 protostars in Orion. 


Manoj et al. 


(120131) noted that these sources 


span two orders of magnitude in bolometric luminosity (0.2 L 0 < Lboi < 28 L 0 ), while 
their CO rotation diagrams show that the CO emission can be characterized by two thermal 
components: warm gas with the rotational temperature of T rot ~ 350 K and hot gas with 
T ro t 700 — 900 K. The rotational temperature of ~ 350 K appears to be universal in 
the mid-J range (14 < J < 24) and independent of the bolometric luminosity. Additional 
observations of the CO rotational emission were obtained as part of the Herschel key 
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program, “ Water in star forming regions with Herschel ” (WISH) observed 18 embedded 


protostars ( Karska e 


sources ([Green et al 


al 


2013), and “Dust, Ice, and Gas In Time” (DIGIT) observed 30 


2013 ). These sources also have properties similar to those observed by 


the HOPS program, i.e., all programs found the universal 350 K component in the mid-J 
CO ladder. 

An LVG analysis indicates that there are two possible explanations for the observed 


CO emission. One is sub-t hermally exci 


low density (< 10 6 cm 3 ) (NeufeL 


d 


2012 


ed gas w ith a high temperature (> 1000 K) and 


Manoj et al 


2013). In this case, one component 


can fit the CO emission over the entire PACS range. The other is thermally excited gas with 
high density (> 10 6 cnr 3 ) and both warm (~ 300 K) and hot (700 — 800 K) temperatures, 
which are shown in the rotational diagram. 


van Kempen et al.l (j2010al) and lYisser et al.l (120121) explored the origin of the mid-J 
CO emission as arising from the UV heated gas along the outflow cavity walls (Photon 
Dominated Region; PDR) and small-scale C-type shocks inside the walls. Furthermore, 


Visser et ah 


(2012) suggested that the PDR contributes more to the FIR CO emission 


as the protostar evolves. However, there is caveats in their model; they combined an 
approximated two-dimensional (2D) PDR model and one-dimensional (ID) shock model, 
and the UV spectrum for gas energetics and chemistry was different in the PDR model. 


Manoj et al 


(2013) argued that the PDR is a minor contributor to the mid-J CO 


emission based on three points. First, it is difficult for the PDR process to produce a 
similar gas temperature regardless of the bolometric luminosity, Tboi- The gas temperature 
in the PDR is roughly proportional to Luy/n, where Luv is the UV luminosity (considered 
to be proportional to L bol ) an d n is the gas density. The density at a specific radius is 
not necessarily correlated with L bo i- As a result, there is no reason to have a similar 
L\jv/n over a large range of Lboi- Second, because of the dilution of UV photons, the 
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PDR cannot produce consistent gas conditions over the large range of radius as required to 
match the CO mid- J emission. Finally, resolved spectra of lower J (< 10) CO emission in 


NGC1333 (Yildiz et ah 


2 01QT) . shows that the contribution of a “broad” spectral component 


is consistent with shocked gas and dominates J— 10-9 emission, implying that the PDR is 


a minor contribution to the CO emission in PACS range. In addition, 


Karska et ah 


(2013) 


also argued that the UV heated gas is a minor component to the mid-J CO emission 
because of a strong correlation between the fluxes of CO J— 14-13 and J— 24-23 and the 
flux of H 2 0 2i2 — loi, which traces shocked gas. 

However, an internal PDR must exist some level in the protostcllar stage. Furthermore 
it may not be negligible because FUV observations toward classical T Tauri stars find 
that t 


(e-g- 


lese stars emi 

the 

JV photons at a 

Herczes et al. 

2002; 

Yang et al. 

2012) 


2012). These UV photons can affect the physical 


and chemical properties of 


(van Kempen et ah 


regard, 


Yildiz et ah 


2(H pa: 


exposed gas with i n the outflow ca vity and along the walls 


Visser et ah 


2012 


Lee et ah 


2014b, hereafter Paper I). In this 


(2012) used spectrally resolved observations of 13 CO J— 6-5 finding that 
the narrow emission lines (An < 2kms -1 ) towards NGC1333-IRAS4A are consistent with 
emission from UV heated outflow cavity walls, which encapsulate the broad outflow lines 
(An > 10 kms -1 ). They find that the mass of the UV heated gas is at least comparable to 
the mass of outflowing gas. In addition, UV photons produced from accreting protostars are 


required to explain t 
with Herschel/HlFl 


the young star (Kristensen et al 


re emission from ionized hydrides (CH + and OH + ) that are detected 


de Graauw et al 


2010) and are inferred to emit from within 100 AU of 


2013). Therefore, a quantitative test with a self-consistent 


PDR model is require to study the importance of PDR for the mid-J CO transitions in the 
embedded protostellar objects. 


We developed a 2D PDR model that self-consistently calculates the gas energetics 
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and chemistry for a given UV spectral type (Paper I). This model was applied to HH46, 
assuming the Draine interstellar radiation held and a T e g = 1.5 x 10 4 K black body radiation 
held (hereafter BB1.5). We find that the PDR reproduces the observed mid- J CO emission, 


while the T e g = 10 4 K black boc 
temperature in agreement with 


y radiation 


Visser et ah 


rereafter BB1.0) results in a lower rotational 


20 121 ). According to this model, the mid-J 


CO emission is radiated from the thermally excited dense CO gas (n > 10 6 cm 3 ) with a 
gas temperature higher than ~ 300 K and the CO abundance above 10~ 5 . 


In this paper, we apply our self-consistent PDR model to a large number of selected 
embedded protostars and test whether the UV heated cavity walls can reproduce the 
universal rotational temperature and huxes observed in the mid-J CO ladder. We present 
properties of our sources in Section [2j and the PDR model and adopted physical parameters 
are described in Section [3j We present our modeling results in Section [4] and discuss 
the effect of physical parameters in Section [51 Finally, we summarize our conclusions in 
Section [6] 


Sources 


In order to test whether the universal rotational temperature can be produced by a 


PDR along the outflow cavity walls, we have applied our PDR model to t 


re sources that 


have pre-existing determinat ions of the density structure in the literature fjjprgensen et al 


2002 


Kristensen et al. 


NGC7129-FIRS2 (Crimicr ct al. 


2012b In addi 


20101 ; 


ion, one intermediate mass embedded protostar, 


Fich et al. 


20101 ) has been modeled as an example of 


a high luminosity source. 


The selected sources (0.8 L Q < Lboi < 500 L e ) are listed in Tabic Q] and are plotted 
in the domain of bolometric luminosity versus the density at 1000 AU (n 10 ooAu) (Figure 
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[[]) . Class I sources (shown as circles) generally have lower values of rtioooAu than Class 0 
sources (squares), as shown in Figure CD We classify sources as “compact” and “extended” 
depending on the distribution of the CO <7=14-13 or CO <7=16-15 emission (in Figured] 
“compact” and “extended” sourc es are marked with open and filled symbols, respectively) 


following the definition of 


Karska et al. 


(20131). 


The observed and synthesized CO fluxes are re prese ntec 


molecules emitting in a given <7 level as follows (IKarska et al 


as the total number of CO 


2o3; 


Green et al. 


20131), 


A/”obs (J) — 


AnD 2 Fj 

hvjAj 


( 1 ) 


where Fj and uj denote the line flux and the frequency of the CO rotational transition from 
J to <7—1, D is the distance to the source, Aj is the Einstein coefficient, and h is Planck’s 
constant. 


Rotational diagrams for our sources are plotted in Figures [2H5] The high-,7 (J > 24) 
CO transitions, which produce a high rotational temperature of ~700-900 K, were detected 
toward most of our sources except NGC1333-IRAS2A and TMC1A. This hot component 
is generally in terp reted as emitting from shocked gas as opposed to the UV heated cavity 


walls (IVisser et al 


2012 


, Paper I). To determine the potential emission from the PDR we 


must first remove the contribution of the hot component from the rnid-J CO emission. 


We calculate two rotational temperatures from the observed rnid-J CO lines. The 
observed rotational temperature of the Warm component, Tw(OBS) is linear-fitted from the 
total observed fluxes, while the Corrected rotational temperature of the Warm component, 
T^(OBS) is derived after subtracting the contribution by the hot component from the total 
rnid-J CO fluxes. For this subtraction, we calculate the mid-<7 fluxes emerging from the 
hot component using the rotational temperature of the Hot component derived from the 
observed high-<7 CO fluxes at J > 24 (T H (OBS)). T H (OBS), T^(OBS), and T w (OBS) for 
each source are listed in Tableland plotted as red, green, and blue color lines, respectively, 




















in Figures EH5J 


We classify sources as TYPE H, P, and S. TYPE H sources, which represent 
half of our sample, are contaminated significantly by the “HOT” component, so 
T^(OBS) < T w (OBS) — 3<j\y(OBS) where ctw(OBS) is a linear fit error of T w (OBS). 
TYPE P (“PURE”) sources, which cover a quarter of our sources, are not contaminated 
by the hot component, so T^(OBS) > Tw(OBS) — 3cxw(OBS). Finally, TYPE S are best 
fit with a “SINGLE” temperature regardless whether it is hot or warm. If the UV heated 
outflow cavity wall reproduces T^-(OBS) as well as the corrected fluxes for the TYPE H 
sources, the hot component is important to produce the universal T rot of 350 K. However, 
for the TYPE P sources, the UV heated outflow cavity walls can be tested directly for the 
universal 350 K rotational temperature. The type of each source is described in Table [2] 
and inside each panel in Figs. EHHl This classification is not correlated with L^ 0 \, Uioooau, 
or evolutionary stage as shown in Figure CD 


The mid-J CO fluxes, corrected for the hot component, still exhibit a nearly constant 
rotational temperature with T^(OBS) ~250 K, as marked with the solid line in Figure El 
(left panel). The right panel of Figure [6] shows the relative contribution of the warm 
component to the total mid-J CO fluxes; more than half of the flux for a given mid-J 
transition are emitted by the warm gas component at J <21. The CO number in J= 14 


from only the warm component fVw(L 
These correlations have been shown in 


) is correlated with nioooAu and Lboi (Figure [7]). 


Karska et al. 


(120131 ). when 1V W (14) was calculated 


from the total flux of J= 14-13 emitted by both warm and hot components. 
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3. Model 

3.1. Density distribution 

We assume that the density in the envelope has a power law distribution of a 
spherically symmetric sphere, excluding the outflow cavity. For our study, the envelope 
ensit y st ructure of each source is dete rmined from previous efforts in the literature 


(Jprgensen et al. 

2002; 

Kristensen et al. 

2012), 


DUSTY (Ivezic & Elitzur 


2012), using the ID radiative transfer program 


19971) . Within this framework the outflo w cavity is carved out 


using the function given below in the Cartesian coordinate system flBruderer et ah 


2009 


z = S 0 x (x 2 + y 2 ) 

1 


( 2 ) 


x (* 2 + y 2 ) 


10 4 AU tan 2 (ct/2) 

where z is the outflow axis and a is the full opening angle at z = 10 4 AU. For the 
ensity inside the ou tflow cavity, we adopt the density of shocked gas, n = 6.3xlO 3 cnr 3 


(Neufcld ct al. 


20091) . which should be the upper limit for the outflow cavity. 


We introduce a new coordinate axis 5 = z/(x 2 + y 2 ) instead of 6 in the spherical 
coordinate system (r, 6) as shown in Figured] While the 9 coordinate describes a circular 
conical surface, the 5 coordinate provides a circular paraboloid. Both PDR and non-local 
thermal equilibrium line radiative transfer models explore scales ranging from ~ 10 AU to 
~ 10 4 AU, resolving the very narrow regions near the outflow cavity wall surface where the 
warm CO gas exists. As the boundary between the outflow cavity and the envelope (d 0 hi 
Eq. [2]) is a point of the 5 coordinate, the (r, 5) coordinates can simply describe the density 


profile of thin layers near the surface (see Figure 16 in Paper 


n 


<5) coordinates through all the procedures except RADMC-3Eo (see below), which does not 


Therefore, we use the (r, 


1 http: //www.ita.uni-heidelberg.de/~dullemond/software/radmc-3d/ 
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provide the coordinate. 


The opening angle is measured by the modeling of existing molecula r line emission 


maps, for example, of 12 CO rotational transitions (e.g., 


Arce fe Sarg ent 


2006). The 


emission distribution toward some sources suggests that the opening angle increases with 
the protostcllar evolutionar y stag e and spreads out from ~10 deg to 100 deg for Class 


0 and I sources flArce & Sargent 


produce the FIR mid- J CO 
cm' 3 ; Paper I; Visser et al 


2006! ). However, if the UV-heated outflow cavity walls 
ines, they should emerge from inner dense regions (n > 10 6 


2012). These regions are within a few arcseconds and are 


smaller than (or comparable to) the beam sizes of mm/sub-mm wave radio telescopes 
even towards the nearby star forming regions. An additional method to determine the 


opening ang 


Furlan et ah 


e is to fit the spectral energy distribution using dust continuum models (e.g., 


20081) . which are model-dependent. The opening angles derived by the latter 


method are generally smaller than (or similar to) those determined by the former method. 
For example, an opening angle of 30° is derived for TMC1 via both methods, while the 
opening angl e of L1551-IRS5 is 10° and 100° by the SED modeling and the CO map, 


respectively flFurlan et al 


2008 


Arce fe Sargents 


2006). Therefore, it is hard to define “an’" 


opening angle for a source. As a result, we assume the opening angle of 30° for all sources, 
which does not change the ID density profile significantly and fits the FIR mid-J CO lines 
reasonably well, compared to other values. In addition, this opening angle produces one of 
the highest CO fluxes for a given UV luminosity, and thus we can test the contribution of 
the UV heated outflow cavity wall to the mid-J CO emission. The effect of the opening 
angle will be discussed in Section 15.11 
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3.2. PDR model 


We have developed a self-consistent PDR model (Paper I). Onr PDR model consists of 
four parts: the calculation of dust temperature, radiative transfer of UV photons, chemistry, 
and gas energetics. The dust temperature T dust is calculated with the dust continuum 
radiative code, RADMC-3D, adopting the dust opacity for the ave rage Milky Way dust in 


dense molecular clouds with Ry=5.5 and C/H = 42 ppm in PAHs (iDrainel 12003!) for a given 
density distribution and a given bolometric luminosity, Lboi- 

The FUV radiative transfer is calculated in order to determine the unattenuated FUV 

inction (Ay) following the method of 


strength G 


(12003 1 and 


and averag e visual ex 


Bruderer et a! 


van Zadelhoff et al. 


(2009a). The FUV radiative transfer is calculated for only one 


representative wavelength with photon energy of 9.8 eV, in the middle of the 6 - 13.6 eV 
FUV band. We then measure the FUV strength Go in units of the Habing field (ISRF; 
1.6 x 10 -3 ergs -1 c um 2 ). W e ad opt the same dust properties used for the calculation of the 
dust temperature (IDrainel 2003). 


In the low-mass classical T Tauri stars, accr etion shocks onto the pro t ostar are 


theorized to produce the observed FUV radiation (Calvet & Gullbrin 


1 998 : 


Inglebv et al 


2 01ll) . while for the intermediate mass Herbig Ae/Be stars, the central star itself can be 


also a FUV radiation source (e.g., 


Meeus et al. 


2012b Bow shocks or small scale-shocks 


inside the cavity or along 

die cavitv wa 

can 

(Neufeld & Dalgarno 

1989; 

Lefloch et al. 

2005) 


so produ ce ad ditional local UV photons 


van Kempen et al 


(2009a) argued that the 


outflow cavity wall could be illuminated by the FUV radiation of ~600 ISRF from the 
shocks in HH 46. In add ition, some sources, such as Elias 29, GSS30-IRS1 (ILiseau et al. 

~ 1 


19991 ). and RCrA IRS5 (iLindberg et al. 


20141 ) are externally illuminated by nearby bright 


stars. 


In our tests, we assume that the only FUV source is accretion onto the protostar. The 
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FUV spectrum affects the photoelectric heating rate of PAHs and small grai ns flSpaans et ah 


1994) as well as photodissociation (and photoionization) of species (van Dishoeck et ah 


2 006 1). However, because we cannot observe the FUV spectrum directly from the central 
protostar, we assume that it is similar to that of a black body radiation of ~15,000 K 


(BB 


.5), which represents the FUV continuum of TW Hya ( Herczeg et al 


2002 


Yang et ah 


2012) and fitted FIR mid- J CO fluxes of HH46 better than a black body radiation of 10,000 


K (BB1.0, Paper I). 

FUV observations toward classical T Tauri stars find the UV luminosity integrated 
from 1250 A to 1750 A (L“) is related with the accretion luminosity L , ]rr as 


log 10 £uv = 0-836 x log 10 L acc — 1.67 with an accuracy of 0.38 dex flYang et ah 
FUV lu minosity integrated f r om 912 A 


2012). As the 


AU Mic (Herczeg et al. 


2002 


Yang et al. 


o 2050 A is about 2 times LJA for TW Hya and 


2012) and the accretion luminosity dominates the 


bolometric luminosity during the Class 0 and I stages, we adopt a reference FUV luminosity 


T y 

-mjv> 


logio-^uv — 0.836 x log 10 Lboi — 1.37. 


(3) 


In this paper, L^ v is used as the unit of Luv- 


When the spectrum in the FUV range is similar to BB1.5, UV photons can be radiated 


from BB1.5 or from t 


K (Nomura & Milla: 


re bremsstrahlung free-free emission with the temperature of ~30,000 


20051 ). If all Lboi is emitted by either of these mechanisms, the 
blackbody radiation and the free-free emission radiate 28 % and 10 % of L bo i in the FUV 
range, respectively. In our model, it is thus difficult to Luv larger than 0.28 L bo i, which is 
similar to the observed la scatter (5 Ly V ) at the lowest L bo i- Therefore, we assume 5 Ly V 
as the upper limit of Luv in our models. 

In ou r model, the gas-ph ase chemical reaction networ k is bas ed on UMIST2006 


database flWoodall et al. 


2007) modified by 


Bruderer et ah 


fl2009btl . For photoreaction 
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rates, we have adjusted the attenuation factor, 7 , following the method of Rollig et al. 
( 2013 ) and calculated the unattenuated photo reac tion rate with the photodissociation 


and photoionization cross sections provided by 


van Dishoeck et al. 


(120061 ). We follow the 


11IUUC1 UI ll^IUilild 

Caza.ux & Tielens 

HUH u 

( 2002 , 

II IIIuC/J 

2004, 

2010 

1 uuoi diiio v icx piiy Diduiptiun a 

) with the sticking coefficient of 

IIU. UIIUIIIIDU1 |JL1UII II UI 

Hollenbach & McKee 


(1197911 . The neutral gas can deplete onto dust grains and evaporate by thermal and 


non-thermal (photon and cosmic ray) events. We also consider electron attachment to 


grain and cation-grain c 

large 

5 x 10 - 17 s " 1 ( 

Dalgarno 

2006) 


We consider important heating and cooling processes described in Rollig et ah (2007). 


We a djust the photoelectric heating rates of PAHs and small g rains (IWeingartner fe Draine 


2001 

) with the correction factor given by 

Spaans et al. 

(1994) 


(11994 1. We also reduce the H 2 


vibrational heating and cooling rate excited by the FUV photons because only UV photons 
in the range of 912 - 1100 A can pump H 2 . We also calculate the H 2 formation heating, 
gas-grain cooling/heating, and atomic and molecular line cooling (for details see Paper I). 

The chemistry and gas energetics are calculated iteratively. It is very time-consuming 
to calculate the chemistry with the full chemical network. In order to reduce the time, 
a subset of the full chemical network has been adopted for gas energetics. As a check we 
compared the_lD PDR models, with a small network with chemical species described in 


Table 1 of 


Woitke et al. 


(120091 ). to identical models with the full chemical network. We 


find that gas temperatures for the two chemical networks are consistent over the range of 
density and FUV strength relevant to this work. The CO abundances near the surface 
(Ay < 1), however, differ from each other by an order of magnitude. Therefore, the iterative 
calculation of the gas energetics and chemistry use the small network, then the chemistry 
with the full network is calculated with the gas temperature determined using the small 





































network. 


3.3. Line radiative transfer 


We have developed a new line radiative transfer code in general grid (RIG). For details, 
refer to Paper I. The most important strength of RIG is the ability to optimize the grid 
coordinates to a given model. RIG works in any coordinate systems, including Cartesian, 
cylindrical, spherical, and (r, 5) coordinates. As described above, the (r, 6) coordinates are 
optimal to model the envelope with outflow cavity walls, and thus, the grid cell number of 
300 in these coordinates (30 in r and 10 in h) provides adequate spatial resolution. The 


best-fit models with a larger number of grid cells (100 in r and 30 in #) , which is com para 


ole 


Visser et al 


( 2012 ), 


to (in r) or higher than (in 6) the spatial resolutions of the model by 
show similar results to the models with the 300 grid cells. 

Collisi onal rate coefficient s for CO ar e adopted from Le iden Atomic and M olecular 


Database 2 


flSchoier et al. 


Visser et al. 


200511 updated bv lYang et al.l ( 20101) and lNeufeldl ( 201211 . Following 


(2012), we fix the non-thermal Doppler width as 0.8 km s 1 and velocity 


distribution as v(r) = 2 km s 1 \fr^Jr with the inner boundary radi us ry n . Because the 


CO ladders in the PACS wavelength range are generally optically thin (IManoi et al. 
the velocity field does not significantly affect the result. 


2013), 


In order to compare to observations, we have synthesized maps of CO spectra, viewed 
at face-on with 0.1" spatial resolution, using a ray-tracing method; these maps are then 
used to predict the number of CO molecules emitting in the J level with Eq. Q] at a given 
pixel. Most of the mid-J CO line produced by PDRs emits within a depth of ~10 AU 
(0".l at 100 pc) from the surface of the outflow cavity walls, which can be represented in 


2 http://home.strw.lcidenuniv.nl/~moldata/ 
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our model due to the optimized 5 grid (see Section l4~2]h We tested the resolution of 0".05 
and found that the difference in simulated fluxes between 0".l and 0".05 resolutions is less 
than 2% in the rnid-J transitions. An edge-on view reduces the mid-J CO fluxes by up 
to 25% as a result of extinction from the dusty envelope. Most of synthesized mid-J CO 
emission arises within the 9 /7 .4 x 9A4 central spaxel of PACS, (which is a few 1000 AU at 
the distance of our sources). 


4. RESULT 


4.1. Best-fit models 


Rotational diagrams from our best-fit models are plotted as black circles in Figures [2]-[5] 
Rotational temperatures derived from the best-fit models Tw (MODEL) are listed in Table [2] 
and inside each panel of Figures [2H5J Most of our best fit models reproduce the observed 
mid-J CO emission for the sources with |T w (MODEL) — T^(OBS)| < 3<rw(OBS) except 


for CedllO-IRS4, VLA 1623-243, and L1551-IRS5. Cedll0-IRS4 and L1551-IRS5 have 
T\v(OBS) > 400 K, and thus, these sources might be heated mainly by shocks, because the 
UV heated cavity wall cannot produce a rotational temperature above 400 K in our models. 
Although VLA 1623-243 has TV(OBS) of 347 K, shocks could also be the main contributor 
to the mid-J_CO e missio n because this source is known to have prominent outflow emission 


(e-g-, 


Bjerkeli et al. 


2012 ). 


Figure 0 shows the best fit Luv (in unit of Lyv) versu s Lboi of sources. From this 
analysis we find that the extended sources have a higher Luv than the compact sources. 
Sources with the upper limit of Luv (NGC1333 IRAS 4B, Ser SMM4, TMC1, and Elias 
29) are all associated with extended emission and our predictions underproduce the mid-J 
CO fluxes. The UV heated outflow cavity wall generally radiates the compact emission (see 
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below), thus, the extended emission is more likely generated by shocks. 

In addition, the sources where we estimate an upper limit to L uv , except TMC1, 

ir. UTTT?T 12 CH f-IOJ 


(San Jose-Garcia et al. 

2013; 

Yilchz et al. 

2013) 


2013). Furthermore, the HIFI 12 CO J=10-9 flux 


is similar to the flux extrapolated from the mid-J lines (see Figure [TUP . T MC1 also shows 


a str ong outflow detected in the PACS [O I] and [C II] lines (IK arska et al 


20131 : 


Lee et al. 


2014al) . Therefore, in this object the mid-J CO emission could readily be produced by 
shocks as opposed to UV radiation. 

However, for a given UV luminosity, BB1.0 generates a lower rotational temperature, 
but higher CO fluxes than BB1.5 (Paper I). When a UV luminosity of 5 L^ v with 
BB1.0 is used, we find a better fit for both rotation temperature and CO fluxes for 
NGC1333 IRAS 4B, TMC1, and Elias 29, as shown in Figure [TQ1 However, the CO 
J= 14-13 flux in Ser SMM4 is still lower than the observed one by a factor of two. 


More than half of the observed CO J = 14-13 flux is radiated from the 


due extended 


(outflow) emission in Ser SMM4 (Karska et al. 


2013 


Dionatos et al 


20131 ). thus, PDRs 


could reproduce only the compact emission near the protostar detected in the central 
spaxel, if any. 

The intermediate mass embedded class 0 protostar NGC7129-FIRS2 has L bo i= 500 L 0 


(Eiroa et al. 

v± * - 

1998; 

Fuente et al. 

2005) 


an accretion rate of M = 10 5 M 0 yr 1 shows 


up to about 4 M 0 (iPalla fe Stahler 


1991 


hat L acc dominates the contribution to L bo i 


birthline, the emission from the sur: 


contributor to Lboi (IPalla fe Stahlei 


ace o 


1991 


1992 ). For the protostar with M* = 5 M 0 at the 
the protost ar with T e g ~ 10 4 K is the primary 
19931 ). and Luv accounts for 7.5 % of 


1992 


Lboi, which is an order of magnitude higher than the best fit Luv (0.5 Ly V ) of the model 
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with BB1.0. Thus, the best fit Luv may allude that most UV photons radiated from the 
protostar might be blocked by the dense material located in the vicinity of the central 
protostar. 

When the UV spectrum of BB1.5 is adopted, the outflow cavity wall heated by 1 Lj ]V 
can reproduce the observed mid-J CO emission. This might imply that the accretion 
luminosity dominates the contribution to Lboi for NGC 7129-FIRS2 (M > 1CT 5 M© yr -1 ) 
and the relation between L acc and Luv derived from the T-Tauri stars could be expanded 
to even intermediate mass protostars. More detailed study, and likely higher resolution 
observations, are needed to distinguish which interpretation is adequate for the PACS 
observation. 

According to our models, the UV-heated outflow cavity walls could reproduce the 
observed “ compact ” mid-J CO emission with or without the hot shock-heated components 
in the case when the observed rotational temperature is below 400 K. We, however, note 
that our model provides only a possible explanation for the mid-J CO emission, and a 
quantitative contribution of the UV-heated outflow cavity wall should be calculated by 
simultaneous modeling of PDR and shocks. 


4.2. Physical and chemical structure of the UV-heated outflow cavity wall 

Figure [IT] shows the physical and chemical properties of the best-fit model for Ser 
SMM1. The properties are plotted along the horizontal distance, Ai£, from the surface of 
the outflow cavity wall for given z-heights, which are marked with horizontal color lines 
and labels in top left panel of Figure fill Other panels show density n (b), the ratio of 
unattenuated (attenuated) FUV strength to density Gq/u (Gdust/^) (c), dust temperature 
T dust (d), gas temperature T gas (e), and CO abundance X(CO) (f). The regions emitting 
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the majority of flux for CO J= 24-23 (filled circles), J— 14-13 (open squares), and both 
transitions (filled squares) are plotted over the layers in the panels. 


The ratio of the FUV strength to density (Go/n) can be used to parameterize the 
dense gas PDR ( n > 10 6 cm~ 3 ) be cause the physical and chemical properties are similar 
for a given Go/n flKaufman et al.l Il999l . Paper I). More directly, photoelectric heating of 
PAHsjmd s mall grains (oc Go n) and gas-grain collisional cooling (oc n 2 ) determine T gas (see 


Visser et al. 


2012 ). CO is destroyed by photodissociation (oc Gon) and forms by two-body 


reactions (oc n 2 ; dissociative recombination and charge transfer) (see Paper I). A higher 
G 0 /n thus gives a higher T gas , but lower X(CO) near the surface. 


In low-mass star forming regions, the power-law index in the density profile is lower 
than two, and Go/n increases toward the center. Therefore, as the z-height is lowered 
(i.e., colors from purple to red in the panel (a) of Figure fTT]) . T gas increases while X(CO) 
decreases near the surface. However, X(CO) near the surface increases from z=500 AU 
downward (cyan line in the panel (f) of Figure 1111) to equatorial plane (red line). Near log 
Go/n ~ —3, UV photons photodesorb H 2 0 ice into the gas phase preventing all oxygen 
from frozen onto dust grains, and a high T gas (> 300 K) makes the formation rate of CO 
high enough to keep X(CO) high in the inner dense regions (see Figure flTTl. where the FIR 
mid-J CO lines form (Paper I). Both H 2 O photodesorption and the fast CO formation at 
>300 K seem important for the physical and chemical conditions in the embedded phase 
(see Paper I for the detail discussion). 


We find that most of mid-J CO fluxes in the best-fit models are produced within 
specific conditions. The CO J— 24-23 line forms in the central spaxel with log Go/n > —4.5 
and log G 0 > 3 for all our sources (see Figure HU). These regions have a density of 
log n (cm~ 3 )> 6 and depth from the outflow surface wall of A R < a few AU (average visual 
extinction of 0.1 < (Ay) (mag) < 1), where T gas > 300 K and log X(CO) > -5. The CO 
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J= 14-13 line emits from the same gas where CO J— 24-23 forms, but also in the gas with 
T gas ~ 100 K, which is located at a higher (Ay) and a larger distance from the protostar. 

Although some sources show that the CO <7=14-13 emission is radiated from outside 
the central spaxels (see the blue squares in Figure [TT]), most of the CO emission is radiated 
from near the protostar in our models. Therefore, the extended rnid-J CO emission cannot 
be reproduced by the UV heated gas, and is likely associated with shock-heated gas. In this 
work, we used the fluxes extracted over the whole PACS spaxels, and thus, for the extended 
sources, a higher best fit Luv is required to reproduce the total flux. 


5. Discussion 

5.1. Effect of physical parameters 

In order to test the effect of physical parameters, we use the model of LI 157, which 
is a compact source located near the median position of the density and the bolometric 
luminosity plot (Figured])- We explore the effect of UV luminosity, opening angle, and 
power index in the density distribution. We set the standard UV luminosity of protostar as 
2.4 v which is the best-fit value for L1157. 

The effect o/Luv- Figure IT21 shows the effect of UV luminosity in the range of 0.0 
< Luy/Lfy < 100. Luy/Lfy > 10 is unrealistic because Luv cannot exceed 0.28 Lboi 
equivalent to Luv/^uv = 9 f° r BB1.5 (see Section lT2l) . but we test two higher ratios 
of 50 and 100 to cover a high dynamic range of Luv/'U-ioooau • A higher protostellar UV 
luminosity produces a larger number of the CO molecules in a given mid-J level. Thus, T rot 
increases with Luv up to Luv < 10 Ly V , then it is nearly constant with T vot ~300 K. 

The left panel of Figure fl3l shows Go estimated for the surface of the outflow cavity as 
a function of Luv- The density power law index is 1.6 and G o follows the inverse square 
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law of distance from the protostar, thus Go/n increases toward the central source. As 
mentioned in Section 14.21 a dense region (log n > 6) exposed to the UV radiation with 
log G 0 > 3 and log G 0 /n > -4.5 can produce the warm CO gas with T rot ~ 300 K traced by 
the CO mid-J transitions. As Luv increases, the region satisfying this condition expands. 
In addition, the filled circles in the left panel of Figure CDU where the majority of CO 
J=24-23 emission radiates, moves outward. This is because the lowest density satisfying 
the condition decreases with Luv, and there is more mass in the envelope at lower densities. 
Furthermore, following them as Luv increases, the filled circles move horizontally along 
log G 0 ~ 4 down to log n ~ 6.5. When L uv is higher than 5 Ly V (yellow line), the largest 
CO J=24-23 emission is radiated around the 1000 AU. The critical density of J— 24-23 
is log n— 7, and thus, CO in lower densities is not thermalized, resulting in decreases 
J= 24-23 emission in the less dense gas beyond 1000 AU. 

The UV heated cavity walls consist of a lower temperature gas component with 
T gas ~ 100 K as well as the warm gas that we are interested in. Therefore, the CO 
,J— 14-13 transition also traces the cool gas (Paper I). When Luv is low, the majority of 
CO J— 14-13 line emission arises from this cool gas, but the contribution of the warm gas 
to the CO J— 14-13 emission increases with Luv, as shown in the right panel of Figure fl3l 
As a result, T rot increases with Luv- 

However, once Luv is high enough to populate the CO mid-J levels consistent 
with T rot 300 K, then it remains at this rotational temperature even with above this 
fiducial Luv- Exploring this more deeply, the lowest and highest mid-J CO transitions 
(J= 14-13 and J— 24-23) are radiated within a narrow temperature range of 200-400 K 
and 300-600 K, respectively, as shown in Figure [13j This is an intrinsic property of 
a PDR; UV photons heat the gas, but also destroy CO. For example, in the ID PDR 
model (Paper I), much of the J— 24-23 flux is emitted within gas with Ay ~0.3 when 
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log n— 6.5 and log Go= 4 (blue lines in Figure fl4l). If this gas is exposed to a higher UV 
flux of log G 0 = 5 (red lines in Figured!]), then the gas temperature increases, but the CO 
abundance decreases at Ay ~ 0.3. As a result, the center of emission arises from deeper 
layers (Ay ~ 1 ). Here the gas temperature is similar to that of the gas with Ay ~ 0.3 
under log Go= 4. Thus there is some similarity in the excitation conditions even within a 
changing external radiation field. 


The effect of opening angle. As the opening angle increases, more UV photons escape 
and do not interact with the envelope. As a result, at a given UV luminosity, the FUV 
strength along the wall of a larger outflow cavity is decreased, resulting in reduced CO 
emission and T rot in the PACS range (see dotted lines in Figure fl5l). As seen in Figure [I5l 
an increase of Luv is required in order to fit the observed fluxes for a wider opening angle. 
Therefore, our models with an opening angle of 30° require the minimum L\jy. Thus, the 
best-fit model with a low L\jy, assuming fixed opening angle of 30°, can have an improved 
solution with a larger opening angle and a higher L\jy. For L1551-IRS5, a model with the 
upper limit of Luv = 5 L} nr and a . large opening angle of 100 °, which has been derived 


from the CO map (lArce fe Sargent 


20061 ). provides an improved fit to the data. However, 


models with different opening angles in Figure fl5l will change dust continuum image and 
spectral energy distribution. Thus, the density profile should be adjusted when assuming a 
different opening angle, which is not accounted for in our models. 


The effect of the density profile. A variation of the power index changes the amount 
and concentration of dense gas, as well as Gq/u, which could change the overall emission 
profile. In the case of L1157, the power index in the density profile of L1157 has a minor 
impact on the rnid-J line fluxes as shown in Figure fl6l The majority of rnid-J CO emission 
arises from gas near 1000 AU (see the left panel of Figure fT3l) and at this radius changing 
the power index does not alter the emission appreciably. However, in the case where rnid-J 
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emission is generated in gas far at greater distances, then the power index significantly 
affects the result. For example, TMC1 has a density at 1000 AU similar to that of TMC1A, 
and a lower L^ 0 \ than TMC1A only by a factor of 3. The power indexes of TMC1 and 
TMC1A are 1.1 and 1.6, respectively. The best-fit Luv (in the unit of L & ) for TMC1 is 
larger than that for TMC1A at least by an order of magnitude, because the lower power 
index of TMC1 reduces the size of the dense region relevant to the population of mid-J 
CO levels. Therefore, a higher Luv is needed to heat the reduced volume and produce the 
observed mid-J CO emission. 


5.2. PDRs as a candidate mechanism to produce mid-J CO emission 


Manoj et al 


(2013) argued that the PDRs cannot reproduce the constant Tw(OBS) 


because T gas is roughly proportional to Tuv/^ioooau- However, in our modeling, T rot is 
nearly constant as ~ 300 K when log Tuv/^ioooau A —6 as shown in Figure |TTJ In 
addition, the variation of T rot is within ~100 K at log Tuv/Aioooau < —6. Therefore, T rot 
can be considered nearly constant with a scatter of 100 K when the mid-J CO emission 
radiated from the UV heated outflow cavity walls. 


As shown in Figure dH most sources radiate the mid-J CO emission within ~1000 AU. 
Therefore, it is impossible that for an internal PDR to account for any extended CO 
emission. About a half of our sources have compact CO emission (see Table [2]). In addition, 
even extended sources also show that a significant mid-J CO flu x (or emission peak) is 


detected near the protostar (Karska et al. 

2013; Green et al. 

2013). Thus, it is possible that 

some portion of the mid-J CO emission 

jroduced by a PDR, which appears to be the case 

. 

for lower-J 13 CO lines (Yilcliz et al. 2012) 



12 CO J=10-9 HIFI data (San Jose-Garcia et al. 2013: 

Yildiz et al. 2013) is the 
























23 


highest transition where pub 


Ser SMM1 (IKristensen et al 


ished data exist that is spectrosc o pical 


2013). According to 


Yildiz et al. 


y resolved, except for 


(2013), the contribution of 


the narrow velocity component to the total J= 10-9 flux seems to have no correlation 
with the evolutionary stage and varies from zero to 100 % with a median of 42 %. I 11 
addition, the contri bution of the broac 


(J < 10) incr eases 


Manoj et al 


(Ivan Kempen et al. 


velocitv componen 

increases a 

O 

O 

O 

CM 

Yildiz et al. 

2012 , 

2013) 


2013). As mentioned by 


(2013), it is likely that a broad (shocked) component has a larger contribution 


to the mid-J CO emission above that seen for J— 10-9. 


However, a UV-irradiated sho ck might be a contributor to the FIR line emission 


toward some embedded protostars (IGoicoechea et al. 


2012 


Kristensen et al. 


2013 


Lee et al. 


2014a[h A portion of the outflowing gas exists between the protostar and quiescent outflow 
cavity walls and, depending on the presence of dust in the cavity, this gas will be irradiated 
by the accreting star, or perhaps by highly excited H or H 2 locally in the gas itself. Thus, 
although it seems likely that the mid-J emission arises in an outflow, it is possible, or even 
probable, that the UV photons are a key additional contributor, beyond the shock itself, to 
heating outflowing gas. 

CO ladder emission must carry a tremendous amount of information on the physical 
state and evolution of protostars. In gas that is not thermally coupled to dust grains (i.e. in 


PDR or shock heated gas) CO emissions are the main g as coo lant anc 


the energy output from star formation. This work and 


Visser et al. 


thus we are tracing 


(2012) have shown 


that an internal PDR can contribute to mid-J CO emissions, even in the case without a 
contribution from shock heating. However, the broad line width detected for resolved CO 
lines (J > 10) implies that this gas is not quiescent, but rather must be entrained in an 
outflow and subject to shock. Evidence for the p resen ce of UV emission is found in the 


detection of ions such as CH + and OH + (Kristensen et al. 


20131) . which require UV radiation 
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in close proximity to the protost ar to b e created. 
water vapor in shocked gas (e.g., Santangelo et al 


urthermore, the inferred abundance of 


2 0141 ) is below that predicted by shock 


models requiring some mechanism to reduce the anticipated H 2 0 creation in the shock. A 
likely candidate to reduce the water abundance is UV photodissociation generated by the 
protostar or perhaps produced in the shock itself. Thus, we know that both shocks and UV 
photons are present in this environment and detailed modeling of both, encompassing all 
these constraints, will be required to extract the key gas physics that is operating in the gas 
near stars at their birth. 


6. Summary 


We have modeled a UV-heated outflow cavity wall as mechanism to produce the mid- J 
CO line emission detected by Herschel/ PACS towards protostars. We obtain the following 
results: 

• The UV-heated outflow cavity walls can reproduce the observed compact FIR rnid-J CO 
emission, except for the sources with the rotational temperature above 400 K, alone or 
when combined with a hot shocked component. 

• The mid-J (14 < J < 24) CO emission can be radiated from the surface (0.1 < Ay < 1) 
of a dense (n > 10 6 cm~ 3 ) outflow cavity wall with log Gq/h > —4.5 and log Go > 3, where 
X(CO) > 10~ 5 and T gas >300 K. 

• Under the above conditions, the hRO photodesorption and the CO formation rates are 
high enough to keep CO in the warm gas, resulting in the mid-J CO emission. 


Our results could support the result of 


Visser et ah 


(20121) and the possibility that the 


PDR contributes at some level to mid-J CO emission for sources with the UV luminosity 
derived from the T-Tauri stars. 
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Table 1. Source parameters 


source 

D 

Lbo\ 

Tboi 

p a 

r . b 
' in 

r b,c 

1 out 

nioooAu a ’ 


pc 

Lq 

K 


AU 

AU 

cm” 3 

L1448-MM 

232 

8.4 

47 

1.5 

20.7 

1.9(4) 

3.9(6) 

NGC1333-IRAS2A 

235 

35.7 

50 

1.7 

35.9 

1.8(4) 

1.7(6) 

NGC1333-IRAS4A 

235 

9.1 

33 

1.8 

33.5 

3.4(4) 

6.7(6) 

NGC1333-IRAS4B 

235 

4.4 

28 

1.4 

33.5 

2.7(4) 

5.7(6) 

L1527 

140 

1.9 

44 

0.9 

5.4 

6.5(3) 

8.1(5) 

CedllO-IRS4 

125 

0.8 

56 

1.4 

4.1 

5.7(3) 

3.9(5) 

BHR71 

200 

14.8 

44 

1.7 

24.8 

1.2(4) 

1.8(6) 

IRAS 15398 

130 

1.6 

52 

1.4 

6.2 

6.2(3) 

1.6(6) 

VLA 1623-243 

125 

2.6 

35 

1.4 

4.3 

1.0(4) 

7.7(5) 

L483 

200 

10.2 

49 

0.9 

12.5 

1.3(4) 

5.1(5) 

Ser SMM1 

230 

1.9 

26 

1.3 

31.0 

1.6(4) 

4.1(6) 

Ser SMM4 

230 

1.9 

26 

1.0 

6.8 

1.1(4) 

5.4(6) 

Ser SMM3 

230 

5.1 

38 

0.8 

8.9 

1.1(4) 

1.1(6) 

L723 

300 

3.6 

39 

1.2 

8.4 

2.4(4) 

8.0(5) 

B33 

250 

3.3 

36 

1.4 

9.8 

1.2(4) 

1.5(6) 

L1157 

325 

6.5 

39 

1.6 

14.3 

3.1(4) 

2.0(6) 

L1489 

140 

3.8 

200 

1.5 

8.4 

6.7(3) 

1.9(5) 

L1551-IRS5 

140 

24.5 

105 

1.8 

28.9 

2.6(4) 

1.2(6) 

TMR1 

140 

3.8 

133 

1.6 

8.8 

7.9(3) 

2.1(5) 

TMC1A 

140 

2.7 

118 

1.6 

7.7 

6.9(3) 

2.2(5) 
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Table 1—Continued 


source 

D 

pc 

Lbo\ 

1/0 

Rboi 

K 

p a 

r . b 

1 m 

AU 

r Q ut b,c 

AU 

nioooAu a ’' 

cm -3 

TMC1 

140 

0.9 

101 

1.1 

3.7 

6.7(3) 

1.8(5) 

HH46 

450 

27.9 

104 

1.6 

28.5 

2.3(4) 

1.2(6) 

DK Cha 

178 

35.1 

591 

1.6 

12.0 

9.6(3) 

9.2(5) 

GSS30-IRS1 

125 

14.5 

138 

1.6 

16.2 

1.6(4) 

1.7(5) 

Elias 29 

125 

20.1 

386 

1.6 

16.2 

1.6(4) 

8.3(4) 

RN091 

125 

2.6 

340 

1.2 

6.6 

5.9(3) 

3.3(5) 

RCrA-IRS5A 

130 

7.1 

126 

0.8 

10.1 

1.0(4) 

2.8(5) 

NGC7129-FIRS2 d 1260 

500 

- 

1.4 

100.0 

1.8(4) 

1.0(7) 


Note. — Sources above the horizontal line are Class 0, sources below 


are C lass I . Physical parameters 


from 


Jorgensen et ah 


(2002) and 


V, r in , r out , and nioooAu) are adopted 


Kristensen et al. 


( 20121 ). 


a The power law index and the molecular hydrogen number density 
at 1000 AU for the density structure, n (r) = nioooAu x ( ioooau ) ? 

b Inner ( r in ) and outer (r out ) boundary radii. 

c a(b) = a x 10 6 


Crimier et al. 

2010; 

Fich et al. 

o 

T-1 

o 

CM 





















Table 2. Model results 


Source 

Extent a 

Th(OBS) b 

T%{ OBS) c 

Tw(OBS) d 

T w (MODEL) 6 

L\jv f 

source type® Ref. h 

L1448-MM 

E 

771 

224 

324 ± 22 

240 

4.1 

H 

2,3 

NGC1333-IRAS2A 

C 

- 

- 

518 ±170 

195 

0.4 

S 

1 

NGC1333-IRAS4A 

E 

479 

184 

298 ± 19 

184 

5.0 

H 

1 

NGC1333-IRAS4B 

E 

893 

255 

343 ± 20 

205 

5.0 1 

H 

4 

L1527 

E 

634 

248 

299 ± 45 

219 

1.0 

P 

1,2,11 

CedllO-IRS4 

E 

799 

417 

491 ± 61 

197 

1.0 

P 

1 

BHR71 

E 

659 

225 

376 ± 27 

237 

1.5 

H 

1,2 

IRAS 15398 

E 

525 

246 

281 ± 18 

220 

4.7 

P 

1 

VLA 1623-243 

E 

- 

- 

347 ± 13 

216 

1.5 

S 

2 

L483 

C 

719 

289 

355 ± 26 

259 

1.0 

P 

1 

Ser SMM1 

C 

656 

273 

351 ± 31 

264 

2.0 

P 

5 

Ser SMM4 

E 

689 

196 

257 ± 18 

228 

5.0 1 

H 

1,6 

Ser SMM3 

E 

653 

195 

291 ± 18 

244 

3.0 

H 

1,6 

L723 

E 

581 

219 

350 ± 26 

231 

1.5 

H 

1 

B335 

C 

612 

224 

310 ± 18 

207 

2.0 

H 

2 



Table 2—Continued 


Source 

Extent a 

Th(OBS) b 

^w(OBS) 

c Tw(OBS) d T w (MODEL) 6 L uv f source type g Ref. h 

L1157 

C 

801 

266 

360 ± 23 

223 

2.4 

H 

2 

L1489 

C 

739 

282 

373 

± 

30 

312 

1.5 

H 

1,2,11 

L1551-IRS5 

c 

- 

- 

436 

± 

32 

114 

0.1 

S 

2,11 

TMR1 

c 

745 

227 

328 

± 

28 

309 

2.0 

H 

1,2,11 

TMC1A 

E 

- 

- 

445 

± 

47 

298 

1.6 

S 

1,2,11 

TMC1 

E 

775 

290 

362 

± 

29 

307 

5.0 1 

P 

1,2,11 

HH46 

C 

652 

265 

307 

± 

22 

236 

0.6 

P 

1,9 

DK Cha 

C 

1056 

262 

386 

± 

26 

219 

0.3 

H 

2,10 

GSS30-IRS1 

E 

802 

212 

335 

± 

19 

308 

2.0 

H 

2 

Elias 29 

E 

723 

159 

315 

± 

19 

322 

5.0 1 

H 

2 

RN091 

C 

738 

275 

252 

± 

25 

270 

0.7 

P 

1 

RCrA-IRS5A 

C 

893 

241 

293 

± 

17 

267 

3.4 

H 

2,7 


797 250 336 ± 24 267 1.0 H 8 


NGC7129-FIRS2 



a Spatial extent of CO lines. When most of mid-J CO fluxes are detected in the central pixel, we define it as 
a compact source (C), and the other case, as an extended source (E). For sources in ref. 1, we use the extent of 
the CO J=14-13 emission. For sources in ref. 2, the sources with a smaller extent than a point spread function 
at CO J=16-15 is considered as compact. 

b Rotational temperature of the “Hot” component from the observed fluxes above CO <7=24-23. 

c Rotational temperature of the “Warm” component from the observed fluxes between CO 7=14-13 and 7=24- 
23 in condition of removing the contribution of the “Hot” component to the fluxes in these transitions. 

d Rotational temperature of the “Warm” component from the total observed fluxes between CO 7=14-13 and 
7=24-23 


e Rotational temperature from the modeled fluxes between CO 7=14-13 and 7=24-23. 
f The best-fit UV luminosity in unit of L \- }v . (see Eq. [3]) 


g Type of Rotational Diagram for the observed data (see text) 
h l 


Karska et al 


2013), 2: 


Green et a 


(2013). 3: 


Lee et al. 


mi 


van KeniDen et al. 

(2010a 

) io 

van Kempen et al. 

“ ^ ^ WU-L, 

(2010b 

), 11: 

Lee et al. ( 


(.2014), 8: 


fl2014al) 


Herczeg et al 


Fich et al 


(2012), 5: 
(20101). 9: 


'The upper limit of Luv hi our model. The observed mid-7 CO fluxes of these four sources are not fitted even 


with this upper limit of Luv- 
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TYPE H 
TYPE P 

Elios 29 E 


NGC1333-IRAS2A 

DKChoO □ SerSMMI 

HH46 


□ 


GSS30-IRS1 


L1551-IRS5 

BHR71 


NGC1333-IRAS4A 


L483 

□ 

RCrA-IRS5A J ■ 

O I-1 LI 448 —MM 

Ser SMM3 I—I 

■ LI 157 NGC1333-IRAS4B 


__ L723 

(T)) TMR1 _ B335 


LI489 _ 

RN091 

TMC1A 


VLA 1623-243 


Ser SMM4 


IRAS 15398 


LI 527 


TMC1 


Cedi10-IRS4 


J_I_i_I_I_I_I_L_ 


J_I_I_I_I_I_I_L 


10 ' 


10 ' 


10 


n 


(cm 3 ) 


1000 AU 


Fig. 1.— The molecular hydrogen number density at 1000 AU (nioooAu) and the bolometric 
luminosity (L bo i) of the sources. Class I sources (circle) are located upper left of Class 0 
sources (square). NGC7129-FIRS2 is excluded in this plot. Results in Table [2] are also 
plotted. Open (filled) symbols represent the compact (extended) sources. The color of red, 
green, and blue represents the source type of “hT, “P”, and “S”, respectively (see text). 
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L1448-MM 



E„ p (K) 


NGC1333-IRAS2A 



1000 2000 3000 4000 5000 
E up (K) 


NGC1333-IRAS4A 



E up (K) 


NGC1333-IRAS4B 


LI 527 




Ep P (K) E up (K) 


Cedi10-IRS4 



Epp (K) 


Fig. 2.— The CO rotational diagrams for L1448-MM, NGC1333-IRAS2A, NGC1333- 
IRAS4A, NGC1333-IRAS4B, L1527, and CedllO-IRS4, in units of total number of detected 
CO molecules (see Eq. [T]) divided by degeneracy g. The open red diamonds indicate the 
values derived from the Herschelj PACS observations. The red (“Hot” component) and blue 
lines (“Warm” component) are linear fits to the observed fluxes of the high-J (E up > 1700 K) 
and mid-J (550 K < E up < 1700 K) transitions, respectively. The green lines are fitted to 
the mid-J fluxes after subtracting the contribution of the “Hot” component from the total 
fluxes. Dotted lines represent the sum of the red and green lines. The open black circles 
represent the best-fit model to the corrected mid-J CO fluxes, and the purple line represents 
the linear-fit of the best-fit model fluxes. The rotational temperature T rot derived from each 
color line and the source type (see text) are presented in the upper right of the box. 
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BHR71 



E„ p (K) 
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B335 



Epp (K) 


Fig. 3.— The same as Figure [2] except for BHR71, IRAS 15398, VLA 1623-243, L483, Ser 
SMM1, Ser SMM4, Ser SMM3, L723, and B335. 
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Fig. 4.— The same as Figure [2] except for L1157, L1489, L1551-IRS5, TMR1, TMC1A, 
TMC1, HH46, DK Cha, and GSS30-IRS1 
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El ias_29 RN091 




E„ P (K) E up (K) 


RCrA-IRS5A NGC71 29-FIRS2 




E up (K) E up (K) 


Fig. 5.— The same as Figure [2] except for Elias 29, RN091, RCrA-IRS5, and NGC7129- 


FIRS2 
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Fig. 6.— Left: the bolometric luminosity (L bol ) and the rotational temperature T-^(OBS) of 
the sources. Right: the contribution of warm component to the total flux of the level J. The 
symbol is the same as in Figure Q] except that NGC7129 FIRS2 (diamond) is included. In the 
left panel, solid and dotted lines indicates the median of T$(OBS) (when excluded TYPE S 
sources as presented by the blue color) and T w (OBS), respectively. In the right panel, the 
solid and dotted lines represent median, maximum, and minimum values, respectively. 
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Fig. 7.— Correlation of the corrected CO number in J— 14 with density at 1000 AU (left) 
and the bolometric luminosity (right). The Pearson correlation coefficient r is presented on 
the panel, which is larger than r— 0.55, the value in the confidence level of 3 sigma, indicative 
of a tight correlation between the parameters. The symbol is the same as in Figure [6] 
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0 12 3 

R (xIO 4 AU) 


Fig. 8.— Illustration of (r, 5) coordinates. The h 0 indicates the boundary between the 
outflow cavity and the envelope described with the equation [2] The general d; is described 
with the same equation except for a larger opening angle at z= 10000 AU. 
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Fig. 9.— The bolometric luminosity (Lboi) and the best fit UV luminosity Luv (in unit of 
Ljj v ) of the sources. The symbol is the same as in Figure |6] Red, blue, and green indicate 
the class 0, I, and intermediate mass sources, respectively. For some sources (CedllO-IRS4, 
VLA 1623-243, and L1551-IRS5), a PDR model cannot reproduce the observed T rot . The 
sources are represented as black symbols. 
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NGC1 333—IRAS4B 


Ser SMM4 




E„p (K) E„ p (K) 


TMC1 Elias_29 




E UP (K) E up (K) 


Fig. 10.— The same as Figure El but the 10 4 K black body radiation field was used for 


9 (ISan .Tose-Garcia et al. 


2013 : 


TMC1, and 

51ias 29. 

Yildiz et al. 

2013 

) is 


’lias 29. The emission of the HIFI 12 CO J— 10 


2013) is plotted as a filled red diamond, which 


is a similar to the flux extrapolated from the rnid-J lines. 
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Fig. 11.— The results of the best-fit model for Ser SMM1. Density (b), ratio of FUV 
strength to density (c), dust temperature (d), gas temperature (e), and CO abundance (f) 
distributions are plotted along given horizontal cuts in the envelope of Ser SMM1. Each 
color line indicates the physical values for a given z-height (in the unit of AU), which is 
represented with the same color in the panel (a). A R is the horizontal distance from the 
outflow cavity wall surface. The filled circles, open squares, and filled squares on top of the 
lines indicate the grid cells where emissions of J=24-23, 14-13, and both lines are radiated, 
respectively. In the panel (c), solid and dotted lines indicate the ratio of dust attenuated 
and unattenuated (incident) FUV strength to density, respectively. (Figures 11.1-11.28 for 
the other sources are available in the online supplementary material.) 
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Luv 



E up (K) 


Fig. 12.— The effect of UV luminosity in the model of L1157. The Luv of the standard 
model (yellow line) is 2.4 Ly V (= 0.82 L 0 ). Each color represents the UV luminosity scaled 
to L 0V (see Eq. [3]). Color lines indicate the rotational diagrams of models with different 
UV luminosities, and the observed data are plotted with open diamonds. The rotational 
temperatures, T mt shown in the right top of the panel are the values fitted to the mid- J CO 
lines of 550 K < E up < 1700 K. A vertical dashed line indicates the highest (J= 24) levels 
in the mid-J CO lines, which are relevant to this work. Lyv > lOL^y (thin lines) are tested 
to cover a high dynamic range of TuvAtoooau although they are unrealistic for L1157 (see 
text). 
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Fig. 13.— G o at the outflow cavity wall surface for given L uv ’s (left) and the Cumulative 
contribution of different gas temperatures to the fluxes of J= 14-13 (middle) and J= 24-23 
(right). The color lines are the same as in Figure[[2j The filled circles and open squares in the 
left panel represent the conditions where most emissions of J= 24-23 and 14-13 are radiated, 
respectively. The filled circles move along the grey solid arrow as Luv increases. Vertical 
dotted lines indicate the distance from the protostar. The black dashed lines represent 
log Go/n = —4.5 and log Go = 3. In the middle and right panels, vertical and horizontal 
lines indicate T gas = 300 K and the cumulative contribution of 0.5, respectively. 
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Fig. 14.— Gas temperature (solid) and CO abundance (dotted) of ID model with 
log n (cnr 3 )= 6.5 for a given UV strength of log Gq= 4 (blue) and 5 (red). Thick lines 
indicate the region emitting 70% of CO J— 24-23 (see Paper I). 
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Fig. 15.— The effect of opening angle a in the model of L1157. Each color line indicates 
the rotation diagram from the model with a different opening angle a. Dotted lines show 
the dependence of a at a given UV luminosity (T uv = 2.4 L\ ]V ). However, solid lines present 
the best-fitted Luv (on upper right) at a given a. 
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power index 



E up (K) 


Fig. 16.— The same as Figure [15] except for the power index in density for L1157. 
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Fig. 17.— The ratio of the best-fit UV luminosity to density at 1000 AU (LuvAtoooau) 
and the rotational temperature of the warm component in the best fit model Tw(MODEL). 
Tw(MODEL) increases and has a strong correlation with Tuv/^ioooau (r= 0.72 in the con¬ 
fidence level of 5 sigma) up to TuvAhoooau ~ 10~ 6 ^ ien near ly constant with ~ 300 K. The 
grey line and symbols indicate the result of Figure [12] The symbols are the same as Figure [9] 




